Construction and test of a moving boundary model for negative streamer discharges 
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Starting from the minimal model for the electrically interacting particle densities in negative 
streamer discharges, we derive a moving boundary approximation for the ionization fronts. The 
boundary condition on the moving front is found to be of 'kinetic undercooling' type. The boundary 
approximation, whose first results have been published in [Meulenbroek et al, PRL 95, 195004 
(2005)], is then tested against 2-dimensional simulations of the density model. The results suggest 
that our moving boundary approximation adequately represents the essential dynamics of negative 
streamer fronts. 
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Streamers are growing plasma channels extending in 
strong electric fields through large volumes of matter; 
they determine the initial stages of electric breakdown 
equally in technical and natural processes [E S H, 
Negative (anode-directed) streamers, which are the sub- 
ject of this work, can be described on a mesoscopic level 
by a system of reaction advection diffusion equations for 
electron and ion densities coupled to the electric field. 
Numerical solutions [1, i, 0, i, i, [H [ll d of this 
minimal model reveal that the evolution of the streamer 
channel is dominated by a space charge layer forming 
around the tip. This layer enhances the electric field 
in front of the streamer, which leads to rapid growth 
through an efficient impact ionization by field accelerated 
electrons. It furthermore screens the electric field from 
the streamer bulk. Fig. [1] illustrates the development 
of this space charge layer, starting from a smooth initial 
ionization seed. More detailed illustrations of fully devel- 
oped streamer fronts will be presented in Sec. IIVI In the 
fully developed streamer, the width of the space charge 
layer is much smaller than the radius of the streamer 
head; this separation of scales is actually necessary for 
the strong field enhancement ahead of the streamer and 
the field screening from the ionized interior. It suggests a 
moving boundary approximation for the ionization front 
which brings the problem into the form of a Laplacian 
interfacial growth model. Such a model was first formu- 
lated by Lozansky and Firsov [l^ and the concept was 
further detailed in [l^, [l^ ; solutions of such a moving 
boundary approximation were discussed in [H, • 

In the moving boundary approximation presented in 
[13, [13, [H, [III the space charge layer is replaced by 
an infinitesimally thin interface. It is assumed that the 
streamer moves into a non-ionized and electrically neu- 
tral medium, so that outside the steamer the electric po- 
tential if obeys the Laplace equation 
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FIG. 1: Stages of evolution of the space charge layer sur- 
rounding a streamer in a homogeneous background electric 
field of _Eo = 1 pointing downwards. Plotted is the half max- 
imum line of the negative space charge density a — p at each 
instance of time. Similar plots in 3 spatial dimensions with 
radial symmetry were presented and discussed in [HI . IT^ , the 
present simulation is purely 2 dimensional. 



conducting 

(p — const inside the streamer. (2) 

We immediately note that this latter assumption will not 
be essential for our analysis. The interface separating 
these two regions moves with a velocity w„ that depends 
on the local electric field 



The ionized body of the streamer is modeled as ideally 
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Here and below, superscripts ^ indicate the limit value 
on the interface where the interface is approached from 
the outside or the inside, respectively. 

In earlier work the electric potential was taken 



as continuous on the streamer boundary, Lp 



+ — 



, which 



for an ideally conducting streamer interior implies that 
the moving interface is cquipotential, ip'^ — const. How- 
ever, due to the small but finite width of the physical 
space charge layer this assumption is unfounded. Rather 
in moving boundary approximation ip must be discontin- 
uous at the interface. In we modeled the potential 
difference across the interface as 



(4) 



where the length parameter £ accounts for the thickness 
of the ionization front and n is the exterior normal on 
the interface. Derivation and discussion of this boundary 
condition is the subject of the present paper. In the 
context of crystal growth from undercooled melts such a 
boundary condition is known as kinetic undercooling. 

Clearly the model ([T])-(l4]) is intimately related to mov- 
ing boundary models for a variety of different physical 
phenomena like viscous finge ring (see, e.g., [lO, HH) or 
crystal growth (see, e.g., [22, |23]). To derive the model 
in the context of streamer motion, our starting point 
is the minimal streamer model that applies to anode- 
directed discharges in simple gases like nitrogen or argon. 
Cathode-directed discharges or discharges in composite 
gases like air involve additional physical mechanisms 

In Sec. Ill Al of the present paper we describe briefiy 
the minimal streamer model. If diffusion is neglected the 
model allows for planar shock fronts moving with con- 
stant velocity in an externally applied, time independent 
electric field, and some properties of these solutions are 
recalled in Sec. Ill Bl Based on these results, in Sec. lIII Al 
supplemented by the appendix, we present a rigorous 
derivation of the boundary condition ([4]), valid for pla- 
nar fronts in strong electric fields. The relation of our 
model to other moving boundary models is briefly dis- 
cussed in Sec. IIIIBI The crucial question whether the 
model also applies to curved ionization fronts in weaker 
external fields is considered in Sec. IIVI Based on numer- 
ical solutions of the minimal model in two-dimensional 
space we argue that our moving boundary model indeed 
captures the essential physics of fully developed (nega- 
tive) streamer fronts. Our conclusions are summarized 
in Sec. El 



II. COLLECTION OF SOME PREVIOUS 
RESULTS 

A. The minimal streamer model 

The model for negative streamers in simple non- 
attaching gases like nitrogen and argon as used in [7:, 
H, [13, [13, [3] consists of a set of three coupled partial 
differential equations for the electron density cr, the ion 



density p and the electric field E. In dimensionless units, 
it reads 



V-je 
dtp 
V-E 



m «(|E|), 

W «(|E|), 
p-a, E = 



(5) 
(6) 
(7) 



The first two equations are the continuity equations for 
the electrons and the ions while the last is the Coulomb 
equation for the electric field generated by the space 
charge p — tr of electrons and ions, je is the electron 
particle current which we here take as the drift current 
only 



je = -crE. 



(8) 



(For the effect of a diffusive contribution to the current, 
see a recent summary in section 2 of [23| and Sec. IIVI of 
the present paper.) The current of the much heavier ions 
is neglected. |crE| a(|E|) is the generation rate of addi- 
tional electron ion pairs; it is the product of the absolute 
value of the current |crE| times the effective cross section 
a(|E|) which is taken as field dependent; an old and much 
used form for a(|E|) is the Townsend approximation 
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-1/|E| 



(9) 



but our analysis holds for the more general case of 
«(0) 



da{\E\) 

0, ,L:' > 0, 



d|E| 



1 



a(|E|) = 1-0 — 
|E|»i V IE' 



(10) 



The dimensional analysis reducing the physical equa- 
tions to the dimensionless model defined above can be 
found in many previous papers, e.g., in [l,[il,[ii. We 
only note that the intrinsic time and length scales are de- 
fined in terms of the electron mobility and the effective 
ionization cross section and thus are determined by the 
physics on microscopic scales. 



B. Planar ionization fronts 

We here recall essential properties of planar negative 
ionization fronts as derived in [IG,, 2M ■ We consider ion- 
ization fronts propagating in the positive x direction into 
a medium that is completely non-ionized beyond a cer- 
tain point X f (t) 



a = — p for X > Xf{t), 
fj > for X < Xf{t). 



(11) 



Far ahead the front, the electric field is taken to approach 
a constant value: 



E ~Eo Eo > 0, 

a;— *-|-oo 



(12) 
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where x is the unit vector in the a;-direction. For a pla- 
nar front, E evidently can depend only on x and t, and 
Eqs. (O, (fT2| yield a constant field in the non- ionized 
region 



E = — i?o X, X > Xf. 



(13) 



The planar solution of the model takes the form of a 
uniformly translating shock front moving with velocity 



In the comoving coordinate 

^ = x-XfiO) 



(14) 



(15) 



a discontinuity of the electron density is located at ^ = 
0, while the ion density p and the electric field E are 
continuous. 




FIG. 2: Densities and field in a planar ionization front in a 
far field _Eo = 1 as a function of the comoving coordinate 
5 (|15|1 . The front moves to the right with velocity v — Eo- 
Upper panel: electron density a (solid line) and ion density 
p (dashed line); lower panel: electric field (solid line, axis on 
the left) and space charge density p — a (dashed line, axis on 
the right). For Q!(|E|) the Townsend approximation ^ was 
used. 



Fig- m that similarly has appeared in |25j], illustrates 
the spatial profiles in such a uniformly translating front 
for Eq = I. In the non- ionized region ^ > 0, we simply 
have p ^ — a and E = —Eq. In the ionized region 
^ < 0, the propagating electron front creates additional 
electrons and ions as long as a{E) > 0, therefore the 
density of immobile ions p increases monotonically be- 
hind the front. The electrons move such as to screen the 
conducting interior from the applied electric field. They 
form a layer of nonvanishing space charge p — a that sup- 
presses the field E behind the front. 

Analytically the solution in the ionized region ^ < is 



given implicitly by the equations 

Eo 



a{E) 
p{E) 



En - \E\ 



En 



£(0 Piv) V ' 



(16) 
(17) 
(18) 



In the limit i?o 3> 1, these equations will be further evalu- 
ated below. We note that at the shock front, the electron 
density jumps from zero to 



(T = Eoa{Eo), 
and for ^ ^ — oo it approaches the value 



cr( — oo) 



a{ri) drj. 



(19) 



(20) 



Far behind the front, where E is so small that a(|E|) ~ 0, 
the final relaxation of E and of a is exponential in space: 
~ e~^^, where 



A = 



cr( — oo) 



1 

— / a{r])d7] 
Eo Jo 



1. 



(21) 



III. THE MOVING BOUNDARY MODEL 
A. Construction 

The results reviewed in the previous section (see also 
Fig. [T]) show that a layer of space charge p — a screens 
the electric field from the streamer interior. For strong 
applied electric field —Eo, the thickness of this layer de- 
fines some small inner scale of the front, while on the 
large outer scale, the streamer will be approximated as 
a sharp interface separating an ionized but electrically 
neutral region inside the streamer from the non-ionized 
charge free region outside the streamer; this substanti- 
ates the assumptions underlying the moving boundary 
model treated in [3 • 

Being a shock front solution of Eq. ([S]), the interface 
always moves with normal velocity 



n- (V^)^ 



(22) 



where n is the unit vector normal to the interface point- 
ing into the exterior region; this equation generalizes 
Eq. We recall that + indicates that the expression 

is evaluated by approaching the interface from outside 
the streamer. 

As mentioned in the introduction, in the context of 
electric breakdown the moving boundary model ((T|) , (HI , 
(|22p has been formulated some time ago by Lozansky and 
Firsov 14]. To complete the model, a boundary condition 



4 



at the interface is needed, and in Ref. fl^ continuity of 
the potential at the interface was postulated: 



if' 



const. 



(23) 



However, Fig. [5] clearly shows that crossing the screening 
layer, (p{^) will increase, which amounts to a jump ip~^ > 
tp~ in the interface model. The size of this jump depends 
on the electric field E"*" at the interface and for a planar 
front is easily determined from Eq. (fT8|) . We note that 
in the framework of the PDE-model, (Sec.|II]), tp^ are to 
be identified as 



ip = (^(-oo), ip+=ip{0), 



f = 



(24) 



(25) 



Since E{^) according to P?|) . (fT5|) is a monotonically de- 
creasing function for ^ < 0, we can integrate by parts 



Ed^ = m 



idE. 



E+ 



£.dE 



(26) 



The last identity holds since either E or ^ vanish on the 
integration boundaries. For a planar front, we have E'^ = 
—Eq, but we here keep E'^ to stress the dependence on 
the field at the front position ^ = 0. 

While E{£^) is known only implicitly as ^ = £,{E) in 
Eq. (fT5)) . the partial integration now allows us to evaluate 
the integral explicitly by substituting Eq. ([T5)l in (pS)) : 



Je 



T]-E+ 



drj. 



(28) 



p{r]) is given in Eq. ^TT} . This result explicitly shows 
that in the interface model the potential is discontinuous 
across the boundary, where the size of the discontinuity 
depends on the electric field right ahead of the ionization 
front. 

Evaluating Eq. (j28p . in the appendix we derive bounds 



showing that tp~^ 



-E~^+ const for large \E^\. We 



here present a simpler argument yielding only the leading 
term. It is based on direct evaluation of Eqs. p?|) . ([18]), 
written as 



\E+\ 



\E^ 



a{f]\E+\)df] 

V — 1 dv 



(5) p{i>E+) V 



(29) 
(30) 



We now take the hmit \E"^\ — > cxD in Eq. (|29p . with 
EjE"^ > fixed. The asymptotic behavior ^ of a(|E|) 
yields 

(31) 

Substituting this result into Eq. ([50)1 we find 

Je(^)/e+ ^ 

yielding a purely exponential front profile 

E{i) = E+e^. (32) 

This result means that the exponential decay of the space 
charge layer ([2T|l . that holds far behind the front for all 

£'+1 > 0, for ^ 1 is actually valid throughout the 
complete front up to <^ = 0. Substituting Eq. p2|) into 
Eq. ([251), we find ^+ - ip- = -E+. 

The more precise argument given in the appendix 
shows that {(p~^ — (p~)/\E^ \ decreases monotonically with 

£'+1 and is bounded as 



< 



< 



\E+\ 
2b 

Je^I 



o 



for |E^ 



where 6 > is some constant. The result 

'ln|£;+| 



ip- = -E+ + b' + 



E+ 



(33) 



oo, 



(34) 



follows. It shows that the first correction to the leading 
behavior ip'^ — tp~ ~ —E^ is just a constant, not a log- 
arithmic term. We thus can choose the gauge of (p as 
ip- +b' ^Oto find (p+ = -E+ + 0(ln \E+\/E+). 

The simplest generalization of this result to curved 
fronts in strong fields, 1, suggests the bound- 

ary condition 



ip' 



n.(V^)+, 



(35) 



replacing the Lozansky-Firsov boundary condition (|23p . 
Boundary condition (|35p together with the Laplace equa- 
tion Aip — ^ and the intcrfacial velocity v = Vtp ([22]) 
define our version of the moving boundary model describ- 
ing the region outside the streamer and the consecutive 
motion of its boundary. 



B. Discussion 

The model formulated here belongs to a class of Lapla- 
cian moving boundary models describing a variety of phe- 
nomena. In particular, it is intimately related to the ex- 
tensively studied models of viscous fingering [2^ |2l[ and 
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solidification in undercooled melts [12, ■ In all these 
models the boundary separates an interior region from 
an exterior region, where the relevant field obeys either 
the Laplace equation or a diffusion equation, and the ve- 
locity of the interface is determined by the gradient of 
this field. 



If we replace the boundary condition p5|) by ((23)) : 
93+ — const., the model becomes equivalent to a sim- 
ple model of viscous fingering where surface tension ef- 
fects are neglected. This "unregularized" model is known 
to exhibit unphysical cusps within finite time [1^, W\ . 
To suppress these cusps, in viscous fingering a bound- 
ary condition involving the curvature of the interface is 
used. The physical mechanism for this boundary condi- 
tion is surface tension. As mentioned in the introduction, 
the kinetic undercooling boundary condition ()35p is used 
in the context of solidification. In that case, however, 
the relevant temperature field obeys the diffusion equa- 
tion. From a purely mathematical point of view, our 
model with specific conditions on the outer boundaries 
far away from the moving interface has been analyzed 
in m, m, [s^l . It has been shown [s^] that with outer 
boundary conditions appropriate for Hele-Shaw cells, the 
kinetic undercooling condition selects the same stable 
Saffman- Taylor finger configuration as curvature regu- 
larization. Furthermore it has been proven [H, [2^ that 
an initially smooth interface stays smooth for some finite 
time. This regularizing property of the boundary condi- 
tion (1551) is also supported by our previous and ongoing 
work [ll,|li,[3l|. 



In applying an interface model to streamer propaga- 
tion, an important difference from the other physical 
systems mentioned above must be noted. For the other 
systems mentioned the moving boundary model directly 
results from the macroscopic physics, irrespective of the 
motion of the boundary: a sharp interface with no inter- 
nal structure a priori is present. In contrast, a streamer 
emerges from a smooth seed of electron density placed 
in a strong electric field, and the screening layer that is 
an essential ingredient of the moving boundary model, 
arises dynamically, with properties determined by the 
electric field and thus coupled to the velocity of the 
boundary. The model therefore does not cover the ini- 
tial "Townsend" or avalanche stage of an electric dis- 
charge [ll| that is also visible in Fig. [1] but can only be 
applied to fully developed negative streamers. Further- 
more, being explicitly derived for planar fronts, the va- 
lidity of the boundary condition for more realistic curved 
streamer fronts has to be tested. This issue is discussed 
in the next section. 



IV. CURVED STREAMER FRONTS IN TWO 
DIMENSIONS 

A. Illustration of numerical results for fully 
developed streamers 

We solve the PDE-model (©-(IT]) in two dimensions, us- 
ing the numerical code described in detail in [ISj . (Pre- 
vious simulation work was in three spatial dimensions 
assuming radial symmetry of the streamer, the results 
are very similar.) In the electron current je, besides the 
drift term — crE, a diffusive contribution is taken into 
account: 

je = -crE - ZJVcr. (36) 

This clearly is adequate physically, and on the technical 
level it smoothes the discontinuous shock front. The price 
to be paid is some ambiguity in defining the position of 
the moving boundary, see below. 

Planar fronts with D > have been analyzed in [l^, 
d, for further discussion and illustrations, we refer 
to section 2 of ^23|. It is found that diffusion enhances 
the front velocity by a term 

= 2V^|S+|a(|S+|), 

that has to be added to the velocity drift = \Ei^\- Fur- 
thermore, diffusion creates a leading edge of the electron 
density in forward direction which decreases exponen- 
tially on scale 

This scale has to be compared to the thickness of the 
screening layer for D = 0: £a ^ l/o^d-^'^l)- Foi' large 
\E^\ and small D both the ratios Wn.n/i'n, drift and £0/^01 
are of order y^D/\E+\ <^ 1. This, by itself, does not im- 
ply that diffusion can be neglected since the term DVcr is 
a singular perturbation of the diffusion free model. How- 
ever, in our numerical solutions the main effect of diffu- 
sion is found to amount to some rescaling of the param- 
eters in the effective moving boundary model, see below. 
This is consistent with the observation that for long wave 
length perturbations of planar fronts, the limit D ^ is 
smooth dl]. 

In our numerical calculations, we take D ~ 0.1. For 
Q!(|E|) the Townsend form ^ is used. We start with 
an electrically neutral, Gaussian shaped ionization seed, 
placed in a constant external electric field E = —Eq x. 
We performed runs for 0.5 < Eq < 2. Since the thick- 
ness of the screening layer decreases with increasing Eq , 
higher fields need better numerical resolution, enhancing 
the numerical cost considerably. In view of the results 
shown below we do not expect qualitative changes for 
E„>2. 

The system ([5])-([7]) is solved numerically with a spatial 
discretization of finite differences in adaptively refined 
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grids and a second-order explicit Runge-Kutta time inte- 
gration, as described in detail in J^j, with the difference 
that the method described there was for three dimen- 
sional streamers with cylindrical symmetry and here we 
adapted it to truly two-dimensional systems. The highest 
spatial resolution in the area around the streamer head 
was Aa; = Ay ~ 1/4 for all simulations. 
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FIG. 3: Zoom into the right half of the symmetric streamer 
head at time step t — 120 of Fig. [1] in a background field 
of _Eo = 1. Shown are level lines of negative space charge 
density a — p (at levels 1/3 and 2/3 of the maximal density) 
and arrows indicating the local electric field at the foot point 
of the arrow. 

For external field Eq = 1, Fig. [1] illustrates the tempo- 
ral evolution of the streamer head. We see that initially, 
a screening layer forms out of an ionization avalanche, 
this process is discussed in detail in The width of 

the layer rapidly reaches some almost constant value that 
depends on Eg . The head develops into a somewhat fiat- 
tened semicircle, with the radius increasing with time. 
This stage of evolution will be addressed as the fully de- 
veloped streamer. Fig.[3]zooms into the streamer head at 
this stage, showing lines of constant charge density p — a 
together with electric field vectors. Evidently screening 
is not complete. A small, essentially constant field exists 
behind the streamer head. This is illustrated in Fig. [3] 
that shows the variation of the electric field and of the 
excess charge along the symmetry line y = 0. This figure 
shows also that the spatial positions of the maxima of |E| 
and of (T — p do not coincide precisely; in fact, the max- 
imum of |E| lies within the diffusive leading edge of the 
front; the small width £d of this diffusive layer replaces 
the jump of the electron density for D — 0. We further- 
more observe that the field behind the front is suppressed 
by about a factor of 1/20 compared to the maximal value 
|E| « 2.1, or equivalently to ^Eq/IO. This screening in- 
creases with increasing background electric field Eq : from 
^Eo/7 for Eo = 0.5 to KiEo/20 for Eq = 2.0. (The max- 
imal field enhancement in these cases is <^ 2 _Bo in the 
fully developed streamer briefiy before branching.) 
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FIG. 4: The same case as in Fig.O plotted is now the negative 
electric field (solid line, axis on the left) and the space charge 
density (dashed line, axis on the right) on the axis of the 
streamer. The plot can be compared to the lower panel in 
Fig. El showing a planar front with vanishing diffusion. 



B. Test of the assumptions of the moving 
boundary model 

The moving boundary model is concerned only with 
the exterior region. Recalling the defining equations 

Aip — 0, in exterior region, (37) 

Vn = -n-E+, (38) 
ip+ = -n-E+, (39) 

we note that all explicit reference to the physics in 
the interior is absent, notwithstanding our derivation in 
Sec. mil Now the first of these equations evidently holds 
as soon as the diffusive leading edge of the electron den- 
sity has a negligible space charge density. Also the second 
relation holds for any smooth shock front {D = 0) of the 
PDE- model. The boundary condition ([5^ . however was 
derived only for planar fronts in strong external fields 
Eo > 1. 

To check whether the moving boundary model ade- 
quately represents the evolution of curved streamer fronts 
for small diffusion and external fields of order unity, we 
first have to choose a precise definition of the interface. 
As illustrated in Figs. [3] and 21 the screening layer is 
fairly thin, but nevertheless it involves the two length 
scales Id and £a and thus shows some intrinsic struc- 
ture. We define the moving boundary to be determined 
by the maximum of |E(x,?/)| along intersections perpen- 
dicular to the boundary. In precise mathematical terms a 
parameter representation (xi,{u) , yi,{u)) of the boundary 
obeys the equation 

= h{u) ■ V\E\\^^^^(^u),y=ytiu), 

where 

n(-u) ^ 

\ c?w ' du J 



7 



is the normal to the boundary at point {xb(u)^yb{u)). 
To motivate this choice we note that the moving bound- 
ary model explicitly refers only to E and not to the 
excess charge distribution. In practice we determine 
{xb{u) , yb(u)) by first searching for the maxima of |E| 
along horizontal or vertical intersections, and we then 
iteratively refine the so determined zero order approxi- 
mation. We always follow the boundary up to the point 
where the local value of jE+j equals Eg. This covers 
all the head of the streamer, where the essential physics 
takes place. Fig. [5] shows the resulting boundary corre- 
sponding to the snapshot of Fig. [31 We observe that the 
direction of E is close to, but does not precisely coincide 
with the normal direction on the interface (except on the 
symmetry axis, of course). The boundary is not equipo- 
tential but a small component of the electric field tangen- 
tial to the boundary drives the electrons toward the tip. 
This effect counteracts the stretching of the space charge 
layer perpendicular to the direction of streamer motion, 
(see Fig. [T]), which in itself would lead to a weakening of 
screening. 




100 



FIG. 5: Again the situation of Figs. |3] and ID Shown is now 
the boundary defined by the maximum of E on normal inter- 
sections (solid line) , the local directions of the surface normal 
(dashed lines) and the local electric fields (solid arrows). 

With the so defined interface we have checked that 
Vn depends linearly on within the numerical preci- 



sion, therefore Eq. (pS]) holds, except for an increase of 
the ratio f„ /|n • E+ 1 , which is an expected effect of dif- 
fusion [H, M HI; this effect can be absorbed into a 
rescaling of time. 

The essential test of the boundary condition ((39| is 
shown in Fig. [6l It shows how ip^ depends on |fi • E+| 
along the boundary for several values of Eq measured at 
times where the streamer is fully developed. For each set 
of data we first determine ip~^ at the maximum of E'^ , i.e. 
at 2/ = 0. This constant ^fi^iE^^^) is subtracted from all 
values ip~^ along the interface. Except for the smallest 
external field Eq = 0.5, the plots in Fig. [S] clearly are 



linear within the scatter of the data. Even for Eq = 0.5 
the curvature is very small. (We note that with increas- 
ing Eq, the width Id of the diffusion layer decreases and 
approaches the limiting spatial resolution of the numer- 
ics ^13] . This explains the increasing scatter of the data 
with increasing Eq). Furthermore, as is illustrated in 
Fig. [7] for Eq = 1.0, for a fixed Eq the slope of the rela- 
tion between (p+ and |n ■ E+| does not depend on time. 
Thus these numerical results can be summarized by the 
relation 



(^+ = MEo,t)-e{Eo)n-E^ 



(40) 



Of course, neither the PDE-model nor the moving bound- 
ary approximation depend on the gauge ipo{EQ,t) which 
thus can be ignored. The prefactor £{Eo) can be ab- 
sorbed into the length scale of the moving boundary 
model, with a compensating change of the time scale to 
preserve Eq. ([55)1 . As mentioned above, this rescaling 
also can absorb the enhancement of w„ due to diffusion. 
As a result, the model ([57)) - ([M)l adequately appears to 
describe also fully developed curved streamer fronts. 
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* 
t 



£(,=0.5 
iJo=1.0 



£(,=2.0 



FIG. 6: ip as a function of iJ+ along the interface for differ- 
ent background fields Eq. A gauge constant 'p'^(E^ax) is sub- 
tracted from (p'^ for each data set. For each electric field, the 
data was extracted at times long enough for the space charge 
layer to be fully developed, but always before the streamer 
branches. These times were t = 500 for iJo = 0.5, f = 120 for 
Eo = 1.0, t = 43.6 for Eq = 1.5 and t = 35 for = 2.0. 



We finally note that the parameter £{Eo) decreases 
with increasing Eq, and it is well conceivable that for 
i?o ^ oo it tends to 1, as predicted by our analysis of 
planar fronts. Furthermore this behavior parallels the 
behavior of the thickness of the screening layer, suggest- 
ing the very plausible assumption that it is this thickness 
which sets the spatial scale of the model also away from 
the limit iSg — > oo. 
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t=80.0 
t=120.0 



1.0 1.2 1.4 1.6 1.8 2.0 2.2 

In-E+I 

FIG. 7: i^"*" as a function of in a background field Eq — 1 
for times t = 80 and 120. The slope is the same. 



V. CONCLUSIONS 

Starting from a PDE-modcl of an anode-directed 
streamer ionization front, we have derived a boundary 
condition valid for a moving boundary model of the 
streamer stage of the discharge. Due to the finite width 
of the space charge layer surrounding the streamer head, 
in a moving boundary approximation the electric poten- 
tial has to be discontinuous across the boundary, and 
the boundary condition (|39[) proposed here accounts for 
this jump in a very simple way. Our analytical deriva- 
tion is restricted to planar fronts in extreme external 
fields Eq ^ oo, but the analysis of numerical solutions 
of the PDE-model shows that the boundary condition 
also applies to (two-dimensional) curved ionization fronts 
in weaker external fields. We conclude that the mov- 
ing boundary model adequately represents the evolution 
of negative streamer fronts. This conclusion can also 
be drawn from studies of periodic arrays of interacting 
streamers that show strong similarities with Saffman- 
Taylor fingers and will be presented elsewhere (3^ . 

As with other moving boundary models in two dimen- 
sions, we now are in a position to use powerful conformal 
mapping techniques to analytically attack questions like 
the stability of streamers against branching. Some first 
results can be found in Refs. [T9|. 

The moving boundary model does not explicitly refer 
to the interior of the streamer and thus leaves open 
questions concerning the role of the residual electric 
field and the resulting currents inside the streamer. 
Analyzing such questions within the framework of the 
minimal PDE model should lead to a more detailed un- 
derstanding of the structure of the space charge layer for 
curved fronts and should clarify the physics underlying 
the phenomenological length parameter £{Eq) occurring 
in Eq. (|40p . This problem, which is important for fully 



understanding the physics of the streamer, is left for 
future work. 
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APPENDIX A: BOUNDS ON ip+ - ip' 

Basic for our discussions are the properties of a (77) 
quoted in Eq. pH)) . valid for 77 > 0: 

< a{i]) < 1, lim 0(77) 1, (Al) 

rj — *oo 



da 
drj 



> 0. 



(A2) 



We furthermore add the physically reasonable condition 
that 

0<b = snJ-^^)<oo (A3) 
exists, so that Q!(77) obeys the bound 

a(v) > 1 - -■ (A4) 
We now rewrite Eq. (fT7|) for p{E) as 

p{E)^i\E+\-\E\)p(^,\E+\Y (A5) 



where 



p{z,\E+\)= a{\E+\iz + {l-z)y)) dy. (A6) 



(A7) 



Eq. P5|) for — tp takes the form 



dx 



\E+\ Jo K^,\E+\y 



The assumption (jAip that a{E) < 1 for all E, leads 
directly to the lower bound 



\E+\ 



> 1. 



(A8) 



We note that a better lower bound can be obtain from 
the fact that since a [E'^ (x(l — y) + y)) increases with 
X, the function p{x,E'^) obeys 

p{x,\E+\)<a[\E+\). 

This leads to the improved lower bound 

\E+\ 



a(|£;+|)' 



(A9) 
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illustrating that weak fields cannot be screened since 
\E+\/a{\E+\) typically diverges ioi\E+\^Q. 

To derive an upper bound valid for large fields, we 
assume \E^\ > b and split the integral in Eq. (|A7|) as 



dx 



Pi^AE+\) 



where 



h 
h 



dx 
dx 

E+\ Pix, \E+\) 



(AlO) 

(All) 
(A12) 



b/\E+\ 



By virtue of Eq. (|A2|) . p{x, \E'^\) increases with x, which 
immediately yields the bound 



h < 



E+ 



Evaluating p(0, |£^+|) with the bound (IA4[) on a{?]) yields 
L n -1 



h < 



\E+\ 



1 - 



\E+ 



(l + ln(|£;+|/6)) 



(A13) 



To evaluate I2 we write 



pix, \E^ 



a{\E+\{x + il-x)y))dy 



ai\E+\iy+{l-y)x))dy 



> I ai\E+\x)dy+ / ai\E+\y)dy, 

> 1 - 7TTT-r + 7TTT-r Inx. 



\E^ 



This result yields 



h < 1 



b/\E+\ f' 



In; 



l-6/|i;+|A/|s+l 1-6(1 -lnx)/|i?+| 



dx, 



< l-b/\E^\ 
1 



l-6(l-ln(6/|i?+|))/|i?+| A/lB+i 
1 



In a; dx. 



l-b/\E^ 



(A14) 



Collecting all the results (and recalling 6/|i?+| < 1), we 
found in this appendix that 



< 



1 



b/\E+\ 



l-^(l + ln(|i?+|/6))' 



|£;+| 1-6/1^^ 

(A15) 

which for large \E^\/b leads to the bound ((33)) given in 
the main text. We note, in particular, that — (p~ does 
not contain a contribution of order ln(|i?+|/6), so that the 
leading (constant), correction to </5+ — (p^ — —E^ can be 
gauged away. 



We finally note that Eqs. imply that 

p{x,\E^\) increases monotonically with \E^\, and thus 
that ((y9+ — Lp~)/\E^\ decreases monotonically. 
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